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A truncated modal expansion approach is used to study the spatiotemporal dy- 
namics of a phase-conjugate resonator as a function of Bragg detuning. The numerical 
results reveal a rich variety of behaviors. Emphasis is given to the spatial distribution 
of optical vortices, their trajectories and their relationship to the beam’s spatial coher- 
ence. The limitations of the model are discussed and experimental results are presented 
for comparison with the model’s predictions and assessment of its soundness. 
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1 Introduction 


Nonlinear dynamics in optical systems (e.g.: passive cavities, wave mixing geometries and 
lasers) have been the subject of intense research during the past few decades. The interest 
stems from the many potential applications of such devices in e.g. optical bistability, phase 
conjugation, image processing, etc. All these applications involve the two dimensional (2- 
D) transverse nature of the system in a fundamental way [1]. Yet, although the nonlinear 
and chaotic behavior of low- dimensional systems described by ordinary differential equations 
(ODE) is fairly well understood [2, 3], problems that include the transverse dimensions and 
are described by systems of partial differential equations (PDE) remain a theoretical chal- 
lenge and do not benefit from well established experimental methodologies [4]. Furthermore, 
their numerical analyses require large computational budgets [5, 6]. 

The system considered in this paper is the phase-conjugate resonator (PCR) sketched 
in fig. 1. It consists of a phase- conjugate mirror (PCM) and a self-imaging cavity made of a 
lens and a planar mirror. PCRs have interesting potential applications in image storage and 
processing [7, 8, 9] as well as in laser systems [10]. The particular geometry of fig. 1 was chosen 
because it leads to simple boundary conditions and because certain aspects of its dynamics 
derived numerically [11] as well as experimentally [12, 13] have already been reported in the 
literature, providing a useful comparative reference frame. In addition, the PCR of fig. 1 
may offer a convenient way of studying experimentally some as yet partially understood 
aspects of more generic dynamical systems. For example, the transverse confinement of the 
resonator is conveniently controlled by the size of two apertures which determine the cavity 
Fresnel number. This is akin to controlling the aspect ratio in hydrodynamics and thus offers 
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a possible way of studying the transition between low aspect ratio dynamics involving few 
modes to larger aspect ratio systems where dislocations, vortices and defects play a dominant 
role [14]. An additional motivation for this work is that experimental and numerical studies of 
many different physical systems are needed to determine possible universal behaviors and test 
theoretical conjectures, e.g. the similar role played by topological defects in different physical 
systems and the striking analogies between laser dynamics and hydrodynamics [15, 16, 17]. 

This paper has two aims. The first is to describe and establish the soundness of a 
somewhat arbitrary yet sensible modal decomposition approach to study the dynamics of 
PCRs with modest Fresnel numbers. This is done by comparing numerical results with 
actual experimental data. The hope is that this comparison, although complicated by the 
difficulty in measuring experimental control parameters exactly and by the problematic effect 
of stochastic noise, will prove sufficiently convincing to confirm the usefulness of the approach. 
The second is to illustrate the relationship existing between the dynamics of the optical 
vortices observed in PCRs [12, 13, 17] and the spatial correlation of the oscillating beam. 

The paper is organized as follow. The PCR’s equations and boundary conditions within 
approximations compatible with experimental conditions are stated in section 2.1. Argu- 
ments justifying the modal expansion approach are given in section 2.2 and the method is 
then used to arrive at a set of modal amplitude equations. In section 3.1, a particular set of 
parameter values is chosen and the system is integrated, using the off-Bragg detuning as a 
control parameter. Different dynamical behaviors are identified from local time series, power 
spectra and phase space portraits. Examples of vortices’ motions and spatial correlation 
maps are given in section 3.2. Experimental results are presented and discussed in section 4 


3 


and a summary is given in section 5. 


2 The model 

2.1 Phase-conjugate resonator 

The PCR is sketched in fig. 1. Its active element is an externally pumped photorefractive 
medium acting as a PCM. Maxwell’s equations for the optical field in the medium lead to 
the following coupled equations: 

{flir + ' v + a + A ^ x 'y ' z ' 1 ) = -Q( x ,yi z i t )M x iyi z i t ) e ~ 3bz , (*) 

- + k 2 ■ V + a + ^ | A^(x, y, z , t) = Q(x, y, z, t)A^{x, y , z, t), (2) 

+ *3 • V + a + A 3 (x , y , z, t ) = -Q(x, y, z, t)A 2 (x, y, z, *), (3) 

^-^■7 + k 4 ■ V + a + -^J A m 4 (x , y, z, f) = Q(x, y, z, t)A*(x, y, z, f (4) 

where, in spirit of the slowly varying envelope approximation, the field in the medium is 
written as 

E(x,y,z,t) = ^/li(x,y,z,f)e j( “ t_fc, ‘ f) + c.c., (5) 

» = 1 

v = c/h is the speed of light in the medium, n is the average refractive index, fc, are unit 
vectors pointing in the direction of the four wavevectors in the medium. The wavenumber 
is ki = ujv where u is the angular frequency of the pump and a is the linear absorption 
coefficient. Vr is the transverse gradient operator. The parameter b is the wavevector 
mismatch (momentum mismatch along the z-axis) measuring the departure from the Bragg 
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condition. This off-Bragg parameter is introduced as a means of varying the amount of phase 
transfer between the waves mixed in the nonlinear medium and is later taken as the control 
parameter. 


The equation for the phase grating Q is derived from the standard band transport model 
of Kukhtarev [18] and reads, in the single grating approximation, 


dQ(x,y,z,t) 

t — ^ \-C\Q(x,y,z,t ) = 7 0 C 2 


' /li(x,y,g,t)A;(x,y,z,t)e jtz + A;(x, y, z, QA 3 (x, y, z, t) 
, h{x,y,z,t ) 


( 6 ) 

where Ci = {Ed + E q +jE 0 )/(E\f + Ep +jE p ) and C 2 = [{E q + Ed)/ Ed}[(Ed + j E 0 )/(Em + 
ED+jE 0 )].r is the Debye relaxation time. Ed, E q and Em are characteristic fields (diffusion, 
limiting space charge and drift field respectively) which depend on the medium properties 
and the grating spatial frequency. E 0 is an external or photovoltaic field, I 0 is the average 
light intensity in the medium and j 0 a coupling constant function of the geometry and the 
electro-optic tensor of the medium. 


The boundary conditions for the wave equation in the medium are determined by the 
cavity geometry. The transverse confinement of the cavity is controlled by two square aper- 
tures Hi and H 2 of sides 2 a pcm and 2a m respectively. The first aperture is located near the 
phase- conjugate mirror in the front focal plane of the lens (fig. 1). The second one is in 
the back focal plane of the lens next to the planar mirror. The transverse confinement is 
measured by the cavity Fresnel number defined as 

T = 2a m a pcm /\f, (7) 


where A is the wavelength in vacuum and / the focal length of the lens. 


\ 
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One cavity round trip, starting at the PCM consists of a truncation by Hx and imaging 
in a 4-f afocal system with a pupil H 2 . Thus, the amplitude A 4 is proportional to the 
truncated convolution of A 3 with a point-spread-function which is the Fourier transform of 
the aperture H 2 . Using transverse coordinates normalized to the aperture size a pcm , the 
boundary condition at the entrance face ( z = 0) of the PCM, which is also the plane of the 
aperture Hi becomes 

A,(x, y, 0, t + 2 L/c) = 7\F 2 rect(x/2)rect(y/2) J J ^ A 3 (x', y', 0, f )sinc [.F(x' + x)] (8) 

sine [;F( y ' + y)] dx'dy', 

where rect(x/2) = lfor— 1<x<1 and = 0 otherwise, sinc(x) — sin(7rx)/:rx and r is the 
amplitude reflection coefficient of the planar mirror. L is the cavity optical length and an 
unimportant phase factor has been omitted in eq. 8. 

The experimental conditions pertaining to section 4 are those of a slow photorefractive 
medium (e.g. BaTiOa). In these conditions, equations 1- 4 can be further simplified by 
neglecting the temporal derivative on the left hand sides. The round trip time delay in the 
boundary condition is also neglected for the same reason. 

2.2 Modal decomposition 

The direct integration of equation 1- 6 requires a large computational budget. Furthermore, 
calculations taking only one transverse dimension into account, althougth showing interesting 
dynamical behaviors (e.g. ref. [11]), fail to reveal the crucial role played by the phase defects 
in mediating the spatiotemporal dynamics of the PCR [12, 13, 17]. 
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Modal decomposition is a powerful tool for studying systems of PDEs [19]. This ap- 
proach allows one to reduce an infinite dimensional problem to a finite (hopefully low) 
dimensional problem described by a system of equations for the modal coefficients. The 
success of the method relies on the assumption that the dynamics is fairly well described 
by the interaction of a small number of active spatial modes. The main difficulty resides in 
the choice of appropriate basis functions for the expansion. They must accurately represent 
the spatial structure of the active modes. However, this structure is generally difficult to 
identify a priori and one often relies on the choice of a sensible basis rather than trying to 
identify the exact optimal basis. 

In optical resonators with relatively small Fresnel numbers, these simplifications appear 
to be justified. The dynamics of these systems is indeed dominated by a few active modes and 
the Gauss-Laguerre or Gauss-Hermite modes of the empty resonator seem to be reasonable 
basis functions to replace the optimal ones [20]. Of course such a practice introduces a degree 
of arbitrariness since a small number of modes of the sensible basis may not be sufficient to 
adequately represent the spatial structure. Nevertheless, the resulting enormous reduction 
of computational budget and the physical insight that can be gained justify further tests of 
the validity of this approach. 

The cavity fields A3 and A4 in the nonlinear medium (0 < z < d) are thus expanded in 
series as 

A 3 (x,y, *,<) = £ hmn(z, t)u„(x, z ; Z a )u* n (y, z; z a ) (9) 

m,n 

A,(x,y,z,t) = J2f^rx(z,t)u m (x,z-,z 0 )u n (y,z;z 0 ) (10) 

771,71 

and the Gauss-Hermite modes are chosen as a sensible set for the decomposition: 
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( 11 ) 


u m (x,z-,z a ) = Q [2 m mlw(z)} 1,2 H m (^^j 

_2 L _2 

exp 27^“^Wr + ^ m + l/ 2 ) tan ^ (*/*.) 1 

w 2 (z) 2R(z) 

where H m are Hermite polynomials. The beam radius is w(z) = u; 0 (l + z 2 lzlY f2 , its radius 
of curvature is R(z) = 2(1 + zl/z 2 ). z 0 = irwlnfX is the Rayleigh distance and k = 2vh/X 
is the wavenumber in the medium. The waist is located at the entrance face of the PCM 
(2 = 0 ). 

The pump beams are Gaussian with a waist w p chosen sizably larger than w 0 and for 
simplicity are assumed to make a very small angle with the z-axis. Thus, the pump beams 
are 

Ai(x,y,z,t) = /i(z,t)*/(x,z)tf(y,x) (12) 

A 2 (x, y, 2 , t) = f 2 (z, t)U*(x , z)U*{y, z ) (13) 

where f/(x, 2 ) = u 0 (x, 2 , 2 p) and z p Trw^TiJX, 

The series in equations 9, 10 are infinite and the problem remains infinite dimensional. 
Because of the cavity transverse confinement however, we expect that only the modes with 
indices smaller than some upper bound will take a significant part in the dynamics. The 
cavity Fresnel number defined by equation 7 is a measure of the amount of transverse spatial 
information that the cavity can accommodate ( T 2 is the space-bandwidth product of the 
cavity). It is thus reasonable to assume that only the modes with indices smaller than T will 
play a dominant role in the wave interaction in the PCM since they are the modes which are 
most likely to survive in the cavity. It must be stressed however that choosing which modes 
to include in the dynamics is somewhat arbitrary. Only a comparison with experiment can 
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provide a a posteriori justification of this choice. The important point is that the series 
of equations 9,10 can often be limited to a small number of terms while still giving results 
which agree at least qualitatively with observations, as will be shown in a later section. 


dz 

dz 

dhmni.Z,t) 


dz 

dz 


&f3mn(Z'>t') = /^(z, 


+ 0 = e ^Mz, (*, t ), 


where 


hmn(-2?0 — JJ dxdyQ^x, y, z, t)u m (z, z)tz n (y, z'jU (x,i)C7 (y,z) 


(14) 


Equations 9- 13 are substituted in 1- 6 and use is made of the fact that the Gauss- 
Hermite modes are solutions of the Helmholtz equation: 

(^ + ^) 

and satisfy the biorthogonality relation [21] 

[ [ dxdyu m {x , z)u n (y, z)u*{x, z)u*(y, z) = S mp 8 nq . 

J J — OO 

A final set of differential equations is then obtained for the modal amplitudes: 


(15) 


(16) 

(17) 

(18) 

(19) 

( 20 ) 


is the overlap integral of the pump and the projection of the grating Q onto the mode (m, n). 


The grating equation 6 is likewise projected onto the cavity modes to obtain 
dh mn (z,t) 


dt 


+ Ci h mn (z,t) = 


( 21 ) 


I) E /«(*. <)<W*. t) + fib, t) E ()}, 

Jfc,/ <:,/ 
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with 


„ i [ f +0 ° ”*(*’ £ Mg: *)“”»(*» z K(y> f)Mfi *)1 2 . d (22 \ 

O klmn (z,t)- JJ _ ^ / 0 (x, y,z,t) 

and the total intensity / 0 (x,t/,z,f) is calculated using equations 9, 10, 12, 13. 

The integration of equations 16- 22 is greatly simplified if the interaction length in the 
PCM is shorter than all the beams’ Rayleigh distances: 

d < z 0 , z p . (23) 

In this case, diffraction can be neglected in the medium and the interacting beams modes, 
normalized to the aperture size a pcm , take the simple form: 

u m (x,z;z 0 ) ~ «m(i) = Q ( 2 ^ 1 ) (24) 

with /3 0 — \f^.cLpcm /w a . The pumps are gaussian with a parameter /3 P — \/2 a pcm / w p- 

The boundary conditions for the modal amplitudes at z = 0 are found from equation 8 
using the modal decompositions 9, 10. They read 

f\mn (0, t) = VT 2 £ / 3W (0, t)I km I , n , (25) 

k,l 

with 

hm = J ^ dxHk{x)u m {x), (26) 

and 

Hk{x) = J dx'ujfc(x')sinc[.F(x' + x)]. (27) 

The algorithm described by Solymar in [22] was used to solve the system of equations 16- 22 
with the boundary conditions 25- 27 . Some results are illustrated in the next section. 
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3 Numerical results 


Phase transfer between interacting beams in the four-wave mixing process appears to be a 
desirable condition for the observation of non trivial dynamics with reasonable values of the 
coupling parameter. For example, the onset of instabilities in ring resonators with injected 
signals has been shown to depend sensitively on this parameter [23]. In a material having 
a real coupling parameter and no external field, this phase transfer can be the result of an 
angular mismatch of the pump beams. For the simulation presented in this section, the 
control parameter was chosen to be the phase mismatch bd (off-Bragg parameter), where b 
is the momentum mismatch along z and d is the length of the nonlinear medium. 

In a recent report [13], it was shown that the dynamics of the PCR strongly depends on 
the cavity Fresnel number. Increasing the Fresnel number relaxes the transverse confinement 
of the system and allows higher order modes to take part in the dynamics. As a result, the 
spatial complexity of the beam increases. In particular, the number of vortices increases and 
their motion becomes more complex, eventually leading to chaotic states and in the limit of 
large Fresnel numbers, to spatiotemporal turbulence. The example shown in the next section 
is for a Fresnel number T — 3.5. This is small enough so that only a few modes are expected 
to take part in the dynamics, yet it is large enough to exhibit a wide range of behaviors. 

Experiments carried out with Fresnel numbers between 3 and 4 reveal a transverse 
pattern in which two pairs of defects of opposite charges repetitively nucleate, move along 
symmetrical trajectories, spending much time in some circular zone centered on the beam and 
disappear at boundaries or collide and annihilate each other [13]. Based on these observations 
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and in order to keep the computational budget to a minimum, it was decided to limit the 
expansions in eqs. 9, 10 to four modes only, namely, (m,n) = (0,0), (1, 1),(2,0), and (0,2) 
with a beam waist w 0 — = 1. 

As already mentioned, there is a certain degree of arbitrariness in this choice. However, 
with a Fresnel number of 3.5, it is unlikely that modes with indices higher than 3 can survive 
in the cavity and play a significant role in the dynamics. Indeed, when a mode of higher 
index was initially added to the series, the simulation showed that its amplitude decays 
and vanishes with the other transients. Another reason for limiting the expansion to these 
particular four modes is that it allows the calculation of the vortices positions to be done 
algebraically rather than numerically, thus saving CPU time. Finally, simulations ran with 
different choices of modes, although giving result which differ in details, revealed the same 
generic range of behaviors. 

The other parameters of the model are chosen as follow: 7 d = —10, ad = 0.15, (3 0 = 
y/2,P P = y/5/10, E d = 1 kV/cm,E q = 5 kV/cm,E M = lOOkV/cm, E 0 = 0 and the pump 
beams intensity coefficients are fi(z = 0) = = d) = 0.5. As an initial condition, it is 

assumed that the four modes are excitied from white noise and have initial modal amplitudes 

Unit = 0 ) = 10 - 4 . 

3.1 Local intensity fluctuations 

The local intensity of the cavity field A a measured at the exit face of the PCM was chosen as 
a variable because it is easily accessible experimentally. The local intensity U{x Q , y 0 , d) was 
calculated at ( x 0 ,y 0 ) = (0.75,0.15) in a region where experiments indicate that the intensity 
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fluctuations have a deep modulation. Figure 2 shows short sections of time series of the 
local intensity fluctuations for different values of the off-Bragg parameter. Time in fig. 2 
is normalized to the relaxation time t' of the photorefractive grating (r' = |(f?A/ + Ed + 
iE 0 )/(E D + E q + i£ 0 )|r). 

After the transients have vanished and a stationary state has been reached, the local 
extrema of these time series were recorded and used to construct the bifurcation diagram of 
fig. 3. For small Bragg mismatch (|W| < 2.4) the motion is periodic. Unstable oscillations 
occur in a wide range of parameter values (3.6 < |W| < 2.4). In most of this range the 
motion appears to be quasiperiodic. Although this is not readily visible at the resolution of 
fig. 3, some regions of the bifurcation diagram seem to be more chaotic than others and there 
exist narrow windows in which the motion is periodic, usually showing a large number of 
subharmonics. For \bd\ > 3.6, the motion is periodic again and eventually, with a sufficiently 
large mismatch, the only stable state is I 4 =0. 

The time series of fig. 2 were chosen to illustrate these various types of motion. The 
corresponding power spectra shown in fig. 4 give a fair indication of the dynamics of the 
system. Additional information can be gained by reconstructing pseudo phase space portraits 
such as these shown in fig. 5 which are plots of the imaginary part of the field A 4 (x 0 ,y 0 ,d) 
versus its real part. 

At bd = —2.3, the motion is clearly periodic with a frequency fj slightly larger than 
1/2 relaxation rate of the grating (fig. 4a) and the phase space portrait (fig. 5a) is a limit 
cycle. At bd = —2.45, a second characteristic frequency f 2 , about 1/10 grating relaxation 
rate, appears and modulate the time series (figs 2b and 4b). For this value of the off-Bragg 
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parameter, the two frequencies are incommensurate and the motion is quasiperiodic with the 
phase space portrait shown in fig. 5b. At bd = -2.55, the low frequency modulation deepens 
(fig. 4c) and more energy shifts toward higher harmonics of fi and linear combinations of 
the two main frequencies (fig. 4c). The phase space portrait still shows some structures 
although it will eventually densely fill a region of phase space since the two frequencies are 
still incommensurate. At bd = -2.75, the time series is more irregular, possibly chaotic 
(fig. 2d). There are still two main frequencies but the spectrum is broadened and shows a 
number of additional peaks indicating the presence of a third frequency about 1/100 grating 
relaxation time (fig. 4d). The phase space portrait of fig. 5d is diffuse and show no apparent 
structure. At bd = —3.18, a phenomenon akin of frequency locking occurs. The two larger 
frequencies pull each other as to become rationed. The very low frequency, which seems to 
have appeared only to allow this locking to occur, disappears. The spectrum of fig. 4e shows 
a ratio ft /f 2 =13. The periodic motion with 13 subharmonics is confirmed by the closed loop 
phase space trajectory of fig. 5e. 

3.2 Vortices and spatial correlation 

The local intensity fluctuations offer a convenient means of comparing model predictions 
with experimental data but this information is not sufficient to fully characterize the spa- 
tiotemporal dynamics. The aim of this section is to illustrate the role played by the phase 
defects (vortices) in mediating the dynamics. 

The vortices appear where the real and imaginary parts of the field amplitude vanish 
simultaneously. For example, fig. 6 shows a snapshot of the beam crosssection I 4 (x,y,d ) for 
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bd = —2.45, which exhibits four dark spots. The corresponding phase contour diagram of 
fig. 7 clearly identifies these spots as two pairs of vortices of opposite charges (±1). In time, 
these vortices nucleates, move around, annihilate each other or disappear at boundaries. This 
is the motion of the vortices which gives rise to the local intensity fluctuations described in 
the previous section. This section gives a more detailed illustration of the spatial aspect of 
the dynamics and attempts to identify the spatial features responsible for the characteristic 
frequencies observed in the local times series. 

Figure 8a to e shows the temporal evolution of the field amplitude A+(x — 0.75, y, d) 
along a line across the beam for the same off- Bragg parameter values used in figs. 2, 4, 5. For 
a small offset (bd = —2.3), all points across the line execute synchronized periodic motions 
with a period a bit smaller than two grating time constants (fig 8a). In 2-D, the phase of 
the wavefront is observed to breathe periodically together with the amplitude but the phase 
gradients produced are not steep enough to trigger the nucleation of defects. 

At bd = —2.45 (fig 8b), the periodic motion is transversely modulated with a period of 
the order of 10 grating time constants. This spatial motion gives rise to the quasiperiodic 
time series of fig. 2b. In 2-D, the phase gradients are locally steep enough to tear the 
wavefront at locations where two pairs of defects of opposite charges nucleate. The four 
vortices then travel across the wavefront, pairs of opposite charge collide and annihilate. 
The higher frequency of the time series seems to correspond to the recovery rate of the local 
wavefront after a vortex has moved through it, while the lower frequency corresponds to the 
full cycle of vortices’ nucleation, motion and annihilation. At bd = —2.55, the motion is 
similar to that just described but with deeper modulation (fig. 8c, which is plotted with a 
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different time scale). 


At bd = —2.75 (fig. 8d), the transverse modulation is irregular, possibly due to the 
occurrence of a third incommensurate frequency. In 2-D, in contrast with the previous case 
where the vortices’ trajectories were confined to fairly well defined area, they now seem to 
visit the whole beam crosssection irregularly. At bd — —3.18 (fig. 8e) the motion is periodic 
again. The transverse modulation has a period of ~ 23 grating time constants, which is 13 
times the period of ~ 9/5 time constant of the fast oscillations measured in fig. 8a. 


The remainder of this section gives a more detailed illustration of the relationship be- 
tween the vortices’ trajectories and the spatiotemporal coherence of the beam. The dynamics 
of the fluctuating beam can be characterized by a correlation index distribution. The corre- 
lation index between two points (xi, j/i) and (x 2 , 3 / 2 ) is defined as the maximum, with respect 
to a temporal shift T, of the correlation of the intensity fluctuations at these two points, i.e. 


C(x 0 ,y 0 ;x,y) = MAX 


' < (hjt)- < h >)(I 2 (t + H- < h >) > 

J< /? > - < h > 2 ]*/ 2 [< Il>-< h > 2 ] 1/2 


VT, (28) 


where /,(*) = I{xi,yi,t);i = 1,2. Stationarity is assumed and < ■ > represents a time aver- 
age. The correlation index is high whenever the intensity fluctuations are highly correlated 
even in the presence of a time delay. Such time delays may be expected if the fluctuations 
are due to features (e.g. vortices) travelling across the beam as waves. 


To illustrate the usefulness of this concept, the correlation index distribution with re- 
spect to the origin, C(0,0;x, y), was calculated for the five off-Bragg parameter values used 
in the previous examples. The spatial correlation diagram thus obtained can then be com- 
pared with the vortices’ trajectories. For this, the positions of the vortices were calculated 
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at regular intervals of time by solving for ReA,(x,y,t) = ImA 4 (x,y,*) =0 and represented 
by dots on a 2-D plot. Note in passing that with the particular, apparently arbitrary choice 
of modes adopted for the expansion of Aj, this could be solved algebraically. 

Figure 9 shows the correlation index distribution for a small value of the off-Bragg 
parameter ( bd = -2.3). No vortices appear in the field but the amplitude and phase of the 
wavefront oscillate periodically as if waves were travelling along an annular area centered on 
the beam. 

At bd = —2.45, the motion is characterized by the nucleation of two pairs of defects 
at two ends of a diagonal. The vortices move along circular trajectories and the members 
of opposite charges of each different pairs collide and annihilate near the two ends of the 
opposite diagonal. This motion is then repeated periodically with alternating directions. 
A sequence of interferogram snapshots illustrating this motion is shown in fig. 10. A plot 
showing the distribution of vortice’s positions in time is shown in fig. 11a. The vortices’ 
trajectories are clearly confined to a narrow annular area in the beam. Figure lib shows the 
corresponding correlation index distribution. Its particular shape can be explained by the 
fact that the intensity fluctuations at some location in the beam are due to the passage of 
a vortex nearby, that no vortex ever appears near the origin (the fixed point for calculating 
C) and that the vortices’ trajectories are confined to a narrow annular region of space. It is 
only near this region that the correlation index is expected to drop sharply. The correlation 
index distribution was also calculated with a fixed point chosen at a location frequently 
visited by a vortex. The distribution looked as a negative of fig. 11a. The correlation index 
is high in the annular area visited by the vortices and low everywheree else. This confirms 
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that the four vortices appearing in the pattern are highly correlated, as is expected since 
they nucleate simultaneously and produce intensity fluctuations near their trajectories which 
are only delayed in time. A qualitatively similar behavior of the vortices was observed for 
bd = —2.55, as illustrated in fig. 12a and b. 

An apparently chaotic time series was observed for bd = —2.75. Figure 13a shows that 
the vortices 5 position distribution is diffuse, indicating that their trajectories do not remain 
confined in narrow regions of space. In fact, the defects can visit about any place in the beam 
except near the origin and this is only because of the peculiar set of modes chosen for the 
decomposition. The corresponding correlation index distribution shown in fig. 13b is low in 
all regions visited by the vortices. Figure 13c shows the correlation index distribution with 
a fixed point (ii, j/i) = (1,0.3) frequently visited by a vortex. It shows that even though 
the vortices trajectories have become more random, the intensity fluctuations in the regions 
they visit remain highly correlated. 

At bd = —3.18, an apparent frequency locking was observed in the local intensity time 
series, leading to a quieter dynamics. For this value of the off-Bragg parameter, the wavefront 
phase and amplitude oscillate periodically but no defects nucleates. The regions in space 
where these fluctuations occur are identified as two concentric rings in the interferogram 
snapshots of fig. 14 and in the corresponding correlation index distribution of fig. 15. 

4 Experimental results 

The experimental apparatus has been described in details in a previous report [13]. The PCM 
is an externally pumped single crystal of BaTi0 3 . The source for the pumps is a single mode 
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Ar ion laser (A = 514nm) optically isolated from the setup by a Faraday rotator. The cavity 
ends with a planar dielectric mirror (R=0.95%) and contains a 16cm focal length lens. Two 
intra-cavity circular apertures of diameter d pcm = 0.79mm and d m =0.38mm give a Fresnel 
number T = 3.7. The off-Bragg parameter was varied by tilting the mirror directing the 
pump A\ toward the PCM and its values were calculated as the product of the momentum 
mismatch along the cavity axis with the interaction length. The change in pump overlap 
due to the tilting was measured and found negligible. 

The aim of the experiments is to show by way of illustrations, that the range of behaviors 
predicted by the model and described in section 3 is indeed observed in the physical system. 
It would be unrealistic to expect a quantitative match between the experimental data and the 
results of numerical simulations if only because it is not possible to control experimentally 
which modes are actually taking part in the dynamics. A qualitative match however, would 
at least indicate that the truncated modal expansion approach gives a fair representation of 
the actual dynamics of the PCR. 

Figure 16 shows an interferogram snapshot of the beam exiting the PCM. It exhibits four 
vortices located at the comers of a quadrangle in an arrangement similar to that predicted 
by the model in fig. 6 and 10. 

The off-Bragg parameter was scanned between two extreme values at which the momen- 
tum mismatch is large enough to reduce the gain below cavity losses and prevent oscillation. 
Local time series were recorded within that range and analyzed using power spectra and 
delayed time phase space portraits. Near the center of the covered range of parameter val- 
ues, a region where the motion was clearly periodic could be identified. On both sides of 
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this region, similar sequences of various dynamical behaviors were observed, up to large mis- 
matches at which the cavity stops oscillating. Before reaching these limits, simple periodic 
oscillations at reduced amplitudes were again observed. The sequence of behaviors just de- 
scribed is exactly what the bifurcation diagram of fig. 3 would predict (note that within the 
approximations of the model, this bifurcation diagram is symmetrical around bd = 0 but 
in the actual setup, this symmetry is broken). By analogy with the model, the origin of 
parameter space ( bd T — 0) was chosen as the midpoint of the central region showing periodic 
oscillations. Other values of the off-Bragg parameter mentioned below are relative to this 
point. Some examples of time series are shown in fig. 17 and the corresponding power spectra 
are displayed in fig. 18. These examples were chosen to illustrate the variety of dynamical 
behaviors exhibited by the PCR and to compare them with the characteristic behaviors 
predicted by the model. 

Figures 17a, 18a show the periodic motion expected for small values of the off-Bragg 
parameter. The corresponding delayed time phase space portrait (fig. 19a) shows broadening 
which may give an indication of the level of stochastic noise (although in our case there is no 
sure way of differentiating between a limit cycle broadened by stochastic noise and a strange 
attractor with a natural width due to the dynamics). 

Figures 17b, 18b obtained with bd r = 0.65 show a chaotic state with an irregular 
time series and a spectrum containing at least two or three incommensurate frequencies, 
including the main frequency of about 0.14Hz, and a large continuum, a signature of chaos. 
The exponential decay of the power spectrum (linear on the semilog scale of fig. 18 may 
be indicative of deterministic chaos). The corresponding phase space portrait of fig. 19b is 
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diffuse and featureless. 


Figures 17c, 18c give a typical example of frequency locking obtained at bd T = 1.29. 
The locking at period three is confirmed by the closed-loop phase space trajectory shown in 
fig. 19c. Note that the width of this limit cycle is about the same as that of fig. 19a and is 
probably due to stochastic noise. 

At bd r = 1.94 (figs. 17d, 18d) the motion is quasiperiodic with a spectrum containing 
only two incommensurate frequencies and their linear combinations. Another quasiperiodic 
motion obtained at bd T = 2.16 is shown in fig. 17e, 18e. Here the lower frequency appears as 
sidebands of the higher one. For larger values of the off-Bragg parameter, e.g. bd T > 2.37, the 
motion is periodic again with the same dominant frequency of ~ 0.14Hz but with decreasing 
amplitude. When bd T > 3, the gain is lower than the cavity losses and no oscillation occurs. 

5 Summary and conclusions 

A model based on a truncated expansion of the field into cavity modes has been tested to 
describe the spatiotemporal dynamics of a phase-conjugate resonator. The Gauss-Hermite 
modes of the empty cavity are chosen as sensible basis functions to represent the fields in 
the nonlinear medium. Although the choice of the modes which take parts in the dynamics 
is arbitrary in this method, justifications are given for truncating the series at a mode index 
of the order of the cavity Fresnel number. Verifications of the validity of this choice are 
later given by comparison with experimental data. Numerical examples, using the off-Bragg 
parameter (momentum mismatch) as a control parameter reveal a rich variety of dynamical 
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behaviors illustrated by the local intensity fluctuations and ranging from simple periodic 
oscillations to quasiperiodic motions and chaotic states. 

The modal decomposition method is particularly well suited to study the spatial aspect 
of the dynamics without requiring prohibitive computational budgets. Numerical results in- 
dicates that the spatiotemporal dynamics of the PCR is mediated by the nucleation of pairs 
of defects of opposite charges in the beam and by their subsequent motion and annihilation. 
Maps of the spatial coherence function of the beam are found to be templates of the corre- 
sponding maps of the vortices’ trajectories, establishing a strong correlation between these 
two quantities. The loss of temporal coherence and the onset of temporal chaos in the local 
intensity fluctuations is likewize correlated to the loss of spatial confinement of the vortices 
trajectories and the loss of spatial coherence. 

Experimental data obtained with a PCR using an externally pumped BaTiC >3 photore- 
fractive PCR reveal a range of dynamical behaviors similar to these shown by the simulation. 
These results confirm and complete the findings of earlier experimental studies of the PCR’s 
dynamics at various Fresnel number [13]. The arbitrariness of choice of the modes to include 
in the dynamics is the main shortcomming of the truncated modal expansion approach used 
in this work. Nevertheless it appears to be a convenient and computationally nondemand- 
ing tool to study and predict qualitatively the spatiotemporal behavior of phase-conjugate 
resonators and perhaps of other nonlinear optical devices as well. 

The authors wish to acknowledge the financial support of NASA/Langley for the re- 
search leading to this report. In particular, we recognize the assistance of Sharon Welch and 
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Figure Captions 


Figure 1: Sketch of the phase-conjugate resonator. PCM, phase-conjugate mirror (externally 
pumped photorefractive crystal); M, planar mirror; Hi, H 2 , apertures controlling the Fresnel 
number of the afocal cavity. 

Figure 2: Sections of time series of the local intensity fluctuations at the output of the PCM 
for different values of the off-Bragg parameter bd (numerical simulation). 


Figure 3: Bifurcation diagram of the local intensity fluctuations at the output of the PCM 
with the off-Bragg parameter as control parameter (numerical simulation). 


Figure 4: Power spectra of the time series of fig. 2 showing a variety of motions: (a)periodic, 
(b,c)quasiperiodic with two incommensurate frequencies, (d)chaotic , and (e)frequency locked 
at period 13. 

Figure 5: Phase space portraits corresponding to the motions illustrated in figs 2, 4 showing 
(a)a limit cycle, (b,c)broadened tori, (d)a diffuse torus, (e)a 13-loop trajectory. 

Figure 6: Snapshot of the field magnitude distribution |Aj| at the output of the PCM 
for bd = —2.45, showing four dark spots identified as four vortices in fig. 7 (numerical 
simulation). 
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Figure 7: Phase contours of the wavefront at the output of the PCM revealing the presence 
of two pairs of vortices of opposite charges. 

Figure 8: Temporal evolution of the field magnitude at the output of the PCM along the 
line x = 0.75 beam waist for the same off-Bragg parameter values as used in fig. 2, 4, 5. 
(a)single frequency periodic oscillation, (b)quasiperiodic oscillations with incommensurate 
transverse modulation, (c)similar to b but with deeper modulation (note the time scale 
change), (d)irregular transverse modulation leading to chaos, (e)frequency locked oscilla- 
tions(period 13) with transverse modulation commensurate with the main frequency. 

Figure 9: Spatial correlation index distribution with respect to the fluctuations at the origin 
for bd = —2.3. No defects nucleates in the beam, but the amplitude and phase breathe 
periodically. 

Figure 10: Interferogram snapshots of the wavefront at the output of the PCM at bd = —2.45 
showing the nucleation, motion and annihilation of two pairs of vortices of opposite charges. 

Figure 11: (a)Distribution of the vortices’ positions at bd = —2.45 showing trajectories 
confined to well defined area. (b)Corresponding spatial correlation index distribution with 
respect to the intensity fluctuation at the origin. 

Figure 12: Same as fig. 11 at bd = —2.55. 
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Figure 13: (a)Distribution of the vortices’ position at bd = —2.75 (chaotic motion) showing 
trajectories randomly visiting large area. (b)Corresponding spatial correlation index distri- 
bution with respect to the intensity fluctuations at the origin, (c) same as b but with respect 
to the intensity fluctuations at a point (x=l, y=0.3) frequently visited by a vortex. 

Figure 14: Interferogram snapshots of the wavefront at the output of the PCM for bd = 
-3.18. No defect nucleates but the phase breathes periodically in two concentric annular 
area. 

Figure 15: Spatial correlation index distribution with respect to the intensity fluctuations 
at the origin for bd = —3.18. 

Figure 16: Interferogram snapshot at the output of the PCM revealing the presence of two 
pairs of vortices of opposite charges, (experimental result) 

Figure 17: Sections of experimental time series of the local intensity fluctuations at the 
output of the PCM for different relative values of the off-Bragg parameter, (experimental 
results) 

Figure 18: Power spectra of the time series of fig. 17 showing a variety of motions: (a)periodic, 
(b)chaotic, (c)frequency locked (period 3), (d,e)quasiperiodic. 
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Figure 19: Delayed time phase space portraits of the local intensity fluctuations of fig. 17, 19 
showing (a)a noisy limit cycle, (b)a featureless portrait, (c)a period three limit cycle and 
(d)a diffuse portrait possibly representing a quasiperiodic motion (experimental results). 
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